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Abstract 

We develop a Petrov-Galerkin stabilization method for multiscale convection-diffusion 
transport systems. Existing stabilization techniques add a limited number of degrees of 
freedom in the form of bubble functions or a modified diffusion, which may not suffi¬ 
cient to stabilize multiscale systems. We seek a local reduced-order model for this kind 
of multiscale transport problems and thus, develop a systematic approach for Ending 
reduced-order approximations of the solution. We start from a Petrov-Galerkin frame¬ 
work using optimal weighting functions. We introduce an auxiliary variable to a mixed 
formulation of the problem. The auxiliary variable stands for the optimal weighting func¬ 
tion. The problem reduces to Ending a test space (a reduced dimensional space for this 
auxiliary variable), which guarantees that the error in the primal variable (representing 
the solution) is close to the projection error of the full solution on the reduced dimensional 
space that approximates the solution. To find the test space, we reformulate some recent 
mixed Generalized Multiscale Finite Element Methods. We introduce snapshots and local 
spectral problems that appropriately define local weight and trial spaces. In particular, we 
use energy minimizing snapshots and local spectral decompositions in the natural norm 
associated with the auxiliary variable. The resulting spectral decomposition adaptively 
identifies and builds the optimal multiscale space to stabilize the system. We discuss 
the stability and its relation to the approximation property of the test space. We design 
online basis functions, which accelerate convergence in the test space, and consequently, 
improve stability. We present several numerical examples and show that one needs a few 
test functions to achieve an error similar to the projection error in the primal variable 
irrespective of the Peclet number. 
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1. Introduction 

Existing techniques for solving multiscale problems usually seek a reduced dimen¬ 
sional approximation for the solution space. Many of these multiscale problems with 
high contrast require stabilization due to the large variations in the medium properties. 
For example, in a multiscale convection-dominated diffusion with a high Peclet number, 
besides finding a reduced order model, one needs to stabilize the system to avoid large 
errors [56]. Stabilization of multiscale methods for convection-diffusion cannot simply use 
a modified diffusion and requires more sophisticated techniques. In this paper, we discuss 
a general framework for stabilization, which combines recent developments in Generalized 
Multiscale Finite Element Method (e.g., [33]) and Discontinuous Petrov-Galerkin method 

(e. g ., mmm)- 

We consider a convection-diffusion equation in the form 

- V • (kVu) + b ■ Wu = f (1) 

with a high Peclet number, where k is a diffusion tensor and b is the velocity vector gU!56]. 
Both fields are characterized by multiscale spatial features. Many solution techniques 
for multiscale problems require a construction of special basis functions on a coarse 
grid (31EHEH HSl I30H331 13711331132113311331137] . These approaches include the Multiscale 
Finite Element Methods (MsFEM) [2] [33] 13711331133] and Variational Multiscale Meth¬ 
ods P, [3] E] [9] [23] I37H3011321133] among others. In MsFEM, local multiscale basis functions 
are constructed for each coarse region. Recently, a general framework, the Generalized 
Multiscale Finite Element Method (GMsFEM), for Ending a reduced approximation was 
proposed (TU] EH EH 113112211331133113311331133] . GMsFEM generates a reduced dimensional 
space on a coarse grid that approximates the solution space by introducing local snapshot 
spaces and appropriate local spectral decompositions. However, a direct application of 
these approaches for singularly-perturbed problems, such as convection-dominated diffu¬ 
sion, faces difficulties due to the poor stability of these schemes. Simplified stabilization 
techniques on a coarse grid are not efficient. Indeed, the modification of the diffusion 
coefficient and similar approaches assumes the use of a few degrees of freedom locally 
to stabilize the problem. These approaches do not suffice for complex problems and one 
needs a systematic method to generate the necessary test spaces. 

We use the discontinuous Petrov-Galerkin (DPG) techniques following |T311231 - 1271 123 . 
[58] to stabilize the system. We start with a stable fine-scale finite element discretization 
that fully resolves all scales of the underlying equation 

Au = f. (2) 
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The system is written in a mixed framework using an auxiliary variable as follows 


(3) 

(4) 


Rw + Au = /, 

A T w = 0. 

The variable w plays the role of a test function and the matrix R is related to the norm 
in which we seek to achieve stability. We assume that the fine-scale system gives w = 0, 
that is, it is discretely stable. In multiscale methods, one approximates the solution using 
a reduced dimensional subspace for u. More precisely, 

u or u ~ Qz u . 

i 

The resulting system also needs a reduced dimensional test space, 

w gs ^ or w « 

i 

The stabilization of ([2]) requires appropriate $ and T. We discuss the design of these 
spaces in the following. 

Within the DPG framework, one can achieve stability by choosing test functions w 
with global support HIES!. However, our goal is to design procedures for constructing 
localized test spaces. In this paper, we design a novel test space which guarantees stability 
for singularly perturbed problems such as convection-dominated diffusion in a multiscale 
media with a high Peclet number. To generate a multiscale space for w, we use the 
recently developed theory for GMsFEM for mixed problems [18]. We start by constructing 
a local snapshot space which approximates the global test functions. These snapshot 
vectors are supported in coarse regions and are constructed solving local adjoint problems 
in neighboring coarse elements. The snapshot spaces are augmented with local bubble 
functions. The dimension of the snapshot space is proportional to the number of fine-grid 
edges (i.e., proportional to the Peclet number). To reduce the dimension of this space, 
making the construction independent of the Peclet number, we propose a set of local 
spectral problems. In these local spectral problems, we use minimum energy snapshot 
vectors ra and perform a local spectral decomposition with respect to the AA T norm. 
Our objective is to find a reduced dimensional approximation, w^, of w such that ||tc— Wn\\ 
is small. We can show that the approximation property of the test space is important 
to achieve the stability (cf. [23] )• We note that the least squares approach [o j8j 3ES GlL] 
can also be used to achieve the stability in the natural norm. Contrary to the traditional 
least squares approach, the proposed approach minimizes the residual with some special 
weights related to the test functions. 

We discuss how to construct online basis [[TT], which uses residual information. Online 
basis functions speed-up convergence at a cost proportional to additional multiscale test 
functions, which are computed by solving local problems. In HI. we developed online 
basis functions for flow equations. One can show that by adding online basis functions, 
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the error reduces by a factor of 1 — A m j n , where A m i n is the smallest eigenvalue for which 
the corresponding eigenvector is not included in the coarse test space. That is, 


.online 


.offline 


w - W "|| < C(1 - A min )|k - < 


where C is independent of the mesh size, physical scales, and material properties’ contrast. 

Thus, if we use all eigenvectors that correspond to asymptotically small eigenvalues in the 
coarse test space, it guarantees that with a few online iterations, we achieve stability. We 
observe this behavior in our numerical simulations. Our construction differs from [13]. In 
this paper, we design different coarse spaces for trial and test. Additionally, the mixed 
formulation we present in ([3])-(|4]) involves higher-order partial-derivative operators than 
standard mixed forms. 

Then, we present several relevant numerical examples of multiscale transport prob¬ 
lems. In particular, we consider heterogeneous velocity fields and a constant diffusion 
such that the resulting Peclet number is high. We consider several types of the velocity 
fields. The first class of velocity fields we use are motivated by [ST] and contain eddies 
and channels. The second class of velocity fields, which are motivated by porous media 
applications, consist of heterogeneous channels (layers). In all examples, we consider how 
the appropriate error (which is based on our stabilization) behaves as we increase the 
number of test functions. We observe that one needs several test functions per coarse 
degree of freedom to achieve an error close to the projection error of the solution of the 
span of the coarse degrees of freedom. Moreover, the number of test functions does not 
change as we increase the Peclet number. By using a few test functions, we can reduce 
the error achieved by standard GMsFEM by several orders of magnitude. 

The paper is organized as follows. In Section 2, we present preliminary results and 
notations, which include the problem setup as well as the coarse and fine mesh descrip¬ 
tions. In Section 3, we describe our proposed procedure. Section 4 contains numerical 
results. Section 5 summarizes our findings and draws conclusions. 

2. Preliminaries 

We consider the following problem 

— V ■ (kVu) + b ■ Va = /, in fI 
u — 0, on dfl 

where k and b are highly heterogeneous multiscale spatial fields with a large ratio maxn(6) / minn(/c). 
The weak formulation of this problem is to find u G V = Hq(Q) such that 


a(u,v) = l(v), Vu G V, 


where 
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We start with a fine-grid (resolved) discretization of the problem and define Uh to be 
the fine-grid finite element solution in the fine-grid space Vf ll Ah and fh are the stiffness 
matrix and the source vector on the fine grid, that is, 

(Ah)ij for (pi, (pj ^ 

(. fh)i = l{4>i) for (pi, (pj e V h 


and u h = YP (pi(uh)i with A h u h = f h . 

We introduce an auxiliary variable (a test variable) and re-write the system in mixed 
form. In particular, we consider the following problem. Find ( Uh°,Wh ) € 14 x Vh such 
that Uh G = Yj <Pi{uh)i and w h = Y <Pi{wh)i solve 

( A h Al A h \ f w h \ _ f f h 

{ ai o ) Ur; _ v o 

Since det(Ah) ^ 0, we have = Uh and Wh = 0. Therefore, these two problems have 
the same solution. Our objective is to find a reduced dimensional coarse approximation 
for Wh, which can guarantee that the corresponding ufY is a good approximation to Uh- 




Figure 1: Illustration of coarse neighborhoods and elements. Red designates a coarse element. Green 
designates two neighboring elements that share a common face (used to construct test functions). Blue 
designates the union of all coarse elements that share a common vertex (used to construct trial functions). 


2.1. Coarse-grid description 

Next, we introduce some notation. We use T H to denote a conforming partition of the 
computational domain D. The set T H is called the coarse grid and the elements of T H 
are called coarse elements. Moreover, H > 0 is the coarse mesh size. We only consider 
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rectangular coarse elements to simplify the discussion and illustrations. The methodology 
presented can be easily extended to coarse elements with more general geometries. Let 
N be the number of nodes in the coarse grid T H , and let {xi | 1 < i < N} be the set of 
nodes in the coarse grid (or coarse nodes for short). For each coarse node Xi, we define a 
coarse neighborhood cu, by 


Ui = U {Kj e T H ; x t G Kj}. (5) 

That is, oji is the union of all coarse elements Kj G T H having the coarse node x t (blue 
region in Figure [Tj) . We use two neighboring elements sharing a common face to construct 
the test functions. An example of this region is depicted in green in Figure [TJ 

3. Generalized Multiscale Finite Element Method for Petrov-Galerkin Ap¬ 
proximations 

In this section, we discuss the construction of the multiscale basis functions for the 
trial space V and the test space W. In particular, we show that one needs a good 
approximation for Wh in order to achieve discrete stability. We start by introducing 
some notation and formulating the multiscale Petrov-Galerkin framework we solve. We 
introduce the snapshot space and then the local spectral decomposition used to construct 
the multiscale basis functions. 

To simplify notations, we let 

Au := —V ■ (kVu) + b ■ Vu, 

and 

A*u := —V • (kVu) — V • (bu). 



(a) Trial snapshot basis (b) Offline trial basis 

Figure 2: A trial snapshot basis and the resulting offline trial basis for a given coarse block. 
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3.1. Construction of the multiscale trial space 

3.1.1. Snapshot space 

We solve a local problem with specifically designed boundary conditions to construct 
the snapshot basis functions. For each coarse neighborhood uji, we define a set of snapshot 
functions such that 

m: ap ) ■= - v • (*v<7 p )+ b ■ v<r = o, 

<t> s ™ P U(. x i) = on dun, 

where is the discrete delta function defined on dcoi with respect to the fine grid. The 
local snapshot space for the trial space is defined by V snap (oui) = span{07 aP }- The 
snapshot functions and multiscale basis functions (offline space) are defined in the union 
of coarse elements that share a common vertex (Figure [l] shows a schematic representation 
of the grid while Figure [2] shows a solution snapshot and a multiscale basis function). We 
use to denote the change of basis matrix from the fine-grid space 14 (tuf) to V snap {uj{). 
Here Vh(coi) is the restriction of 14 in a 

3.1.2. Eigenproblem 

To construct the offline trial space, we solve the following eigenproblem 

(A'£ Mp ) T A% ap v j = A jM2 apVj , 


where 


mz v = *TK“ i >i 

and (Xj, Vj ) is the j-tli eigen-pair. In the above definition, Af l and Mff are the restrictions 
of the fine-scale stiffness matrix Ah and the fine-scale mass matrix Mh in oji. We order the 
eigenvalues in increasing order and we use the first m eigenfunctions as the offline trial 
basis functions. Specifically, we define 4/ = $i v j and 14// = span{y/4,/|l 4 j < m, 1 < 
l < N}, where {xi} is the partition of unity. 

3.2. Construction of the multiscale test space 
3.2.1. Snapshot space 

The snapshot space for the test space consists of three components, and is denoted 
as W sn ap = Wf nap + Wf nap + Wf nap . Next, we will give the constructions for W x snapl Wf nap 
and Wf nap . In each coarse block A4, we define Wf nap (Kk) as 

^Lp(Kk) ■= {r nap 6 14,0 (K k ) I A*(r nap ) = Zij in K k for some 6,/ e 14//} 

where Vh(K k ) is the restriction of 14 in K k and Vh,o{K k ) is the subspace of Vh(K k ) 
containing functions that vanish on dK k . The space Wf nap (K k ) contains functions that are 
solution of the adjoint problem on K k with a source term fij and zero Dirichlet boundary 
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(a) Test snapshot basis for W 3 (b) Offline test basis for W 3 

Figure 3: A test snapshot basis and the resulting offline test basis for a given coarse block of W 3 . 

condition. The space Wf nap , is defined as Wf nap = ® k W} nap (K k ). The space Wj nap is 
considered as the space of multiscale bubble functions. We remark that we obtain perfect 
test functions (with perfect stabilization) if the above local problems are solved on the 
whole domain. 

The second space Wf nap is defined as follows. For each coarse block K k , we define 

W 2 snap {K k ) := {4 snap £ V h (K k ) | A*^ snap ) = 0 in K k and 4 snap is linear on E £ d K k j. 

The space Wf nap is defined by Wf nap = ®kWf n ap( E k) H C* 0 ^). Note that this space is 
similar to the classical multiscale finite element space. 

Finally, we give the definition for Wf nap . For each coarse edge E k , we define K(E k ) 
as the set of all coarse blocks having the edge E k . Then, we find ^ n k ap £ Vh(K(E k )) such 
that 


A ^sna P := _ v . („v^n - V • (b^ ap ) = 0 in each K £ K(E k ), 

^) = S ij for a]1 X j e E k , 

<n mE>)\E> = 0. 

In the above system, 5^ is the discrete delta function defined on E k with respect to the 
fine mesh and is zero on the boundary of E k . Then we define Wf nap {E k ) = span{-0^ ap }, 
and W? nap = @ k Wg nap {E k ). (Figure [Tj illustrates a grid and Figure [ 3 ] illustrates a snapshot 
solution and a multiscale basis function). 

We remark that W snap = 14- 

Lemma 1. For each u £ V 0 ff, there exists a test function </> £ W snap such that 

a(v,<j>) = (u,v)i 2, Vv £ V h . 
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( 6 ) 








Proof: Given u G V q r, we assume that (f> G Vh satisfies ({6]). For each K G T H , we define 
^snap e W'snap satisfying 

a ( v ><f> ( snl P ) = (w,v)p, Vv e V h)0 {K), WK G T h . 

Next, we define 0snl P G W s 2 nap + Vfy 3 nap such that 

<t>£L = 0 on 0K, MK G T h . 


Then, we have 


a(v,0- - «p) = 0, Vt, € V h , 0 (K), VA € r" 

with <f> - 0si\L P - ^inap = 0 oil t?K, VA' 6 T H . Therefore, <p - <t>in ap - etna,, is a solution of 
the adjoint problem with zero Dirichlet boundary condition and zero source term in all 
K G T H . Thus, 0 = 0inip + 0snap G I'fynap- 


3.2.2. Eigenproblem 

Among the three parts of the test space, the dimension of is proportional to the 

number of hne-grid blocks and thus proportional to the Peclet number of the problem. 
Consequently, our objective is to reduce the degrees of freedom associated with W^ nap . 
Both the dimensions of W} nap and W] nap are proportional to the number of coarse grid 
degrees of freedom. We consider two different eigenvalue problems to construct the offline 
test space for W^ nap . 

The first eigenvalue problem for fffy nap (Ify): In this eigenvalue problem, we will 
use the edge values of the snapshot solutions. 




V0 


The eigenvalues go to oo as we refine the fine mesh. 

The second eigenvalue problem for Ify 3 nap (ify): This eigenvalue problem is mo¬ 
tivated by da, where we construct minimum energy snapshot solutions and perform a 
local spectral decomposition using the same norms. More precisely, 



>K(E k ) 


(A t v)(A t ^) 


where 0 = arginin^^^^,^ =^ Ek} {.f K(Ek) (A T rjj) (A T 0)}. In this case, the eigenvalues 
are always smaller than 1. 

We will arrange the eigenvalues of the above spectral problems in increasing order, and 
choose the first L k eigenfunctions as the offline test basis functions. The span of these basis 
functions is denoted as The final test space W„ff is defined by W^ nap 0 W ^ nap © . 
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3.3. Global coupling 

We can use the above trial and test spaces to get a reduced system for the multiscale 
solution. In particular, the multiscale solution is computed by solving 

( ® T A h A T h Q Q T A h E \f w ms \ _ ( Q T f h \ ( . 

\ s T ^e 0 0 )■ (7) 

The columns of 0 consist of the computed multiscale test functions while the columns of 
5 consist of the computed multiscale trial functions. 


3-4- Discussion 

Next, we discuss the approximation properties for the test space and how the defini¬ 
tion of this space affects the discrete stability of the resulting method. To simplify the 
discussion, we introduce some notation. Let N and M be the dimensions for the test and 
trial spaces, respectively. Thus, we can write 

f QjfAhA^Qjv QJ^A^Em \ f u>n,m \ _ ( G>Nfh \ /g\ 

\ 0 /V U N G M J \ 0 J 

For simplicity, we denote by 0^ the snapshot matrix that contains all snapshot vectors 
in the test space and similarly for the trial space E^. Therefore, the following statements 
are true: 


* ^ 00,00 0 . 

• Uoo : m is a projection of Uoo i00 onto Em 

'U'OO ,M nHjy^CX),00' 


• Our objective is to find the smallest possible N and M 0 , such that \\u^ G M — -< 

ll^Sfoo — u™ M \\ for any M, when M > M 0 . 

• The inf-sup condition for our discrete saddle-point problem can be written as 

S q P (0^l4e") 1/2 “ C infsu P ^ T M M E M ) 1/2 . (9) 

The inf-sup condition implies that 


a 


infsup 


= inf 


in 


A T et 




«=“M? 


Mi* 


Next, we take u = A T z. Then, the projection of u onto A t Qm is 


n Ar0 » = A T e N ((A T e N ) T A T e N )- 1 (A T e N ) T A T z 

= A 1 Q N (Q r ^AA T Q N )^ 1 Q^ f AA T z. (10) 
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We define the 0 projection in the B norm to be 


n e<B (z) = 0(e) 1 BO)- 1 ® 1 Bz. 


Thus, 


||n A T 0J »|| 2 = z T AA T 0 N (0 T N AA T e N )- 1 0 T N AA T 0 N (0 T N AA T 0 N )- 1 e T N AA T z 
= z T AA T 0 N {0 T N AA T 0 N y l 0 T N AA T z 
— l|n e ^T(z)||^ T . 


( 11 ) 


Also, 

Thus, 


Cinf sup 

If the inf-sup is satisfied, then we have 


u\\g = IrUat- 

|I1© J v,A j 4t(^)||aAT’ 


inf _ 

=A T z,u=E. M q 


z\\aa t 


( 12 ) 


tVN,M ^oo,oo || T ||^jV,M ^oo,oo II — ||^JV,M ^00,00 | T tl 


PG | 

00,00 1 


= 0 + ||m PG oo,M - U ™\\. 


PG 


From here, we have 


I tVoo,M 
,PG 


PG 


■ PG » -< II ' U PG„ - PG 


+ \\ U oo,M M oo,ooll — \\ ur ^oo,M U oo,oo\ 


(13) 


(14) 


Because N — oo, tr °oo,m = u PG oo,m, we get 


II^oo.mII ^ ||« PG oo,m -« pg 

The discrete inf-sup condition can be shown if for any z (e.g., z = A _t Sm9 ), there 
exists z 0 in the space spanned by 0^ (he., z 0 = Ojqz r ), such that 


\z - Z 0 \\ AA T < ^II^IUat, 


(15) 


for some S < 1. In multiscale methods (in particular, in our works P3 HH]), we 
reduce the error in ||z — £o|Ua t by selecting appropriate multiscale spaces (as those 
used herein). In addition, this procedure can be done adaptively. Thus, by selecting 
a sufficient number of multiscale basis functions, we can reduce the error \\z — 2 o||aa t 
and can achieve the stability sought. We do not have rigorous error estimates, but 
study this problem numerically. We emphasize that we need good approximation 
properties in the test space (as in [24]), which is due to the primal formulation and 
the choice of z in (10). 
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3.5. Online test basis construction (residual-driven correction) 

One can use residual information to construct online basis functions. Online basis 
functions use global information and thus accelerate the convergence. In [HI, we discuss 
the online basis construction for flow equations using a mixed formulation. We use the 
local residual to construct an online basis function locally in each non-overlapping coarse 
grid region a;*. 

The offline solution in the fine-scale test space (woo,m,Uoo,m) £ 14 x 14// satisfies 


f A h A T h A h Sju \ ( Woo ,M \ _ f fh \ 

V -llAl o J \ Uoc,M J V 0 / 

and the multiscale solution (iwjv.Mj wtv,m) G VF 0 // x 14 // satisfies 

f QjfA h A)(Q N Q^A^m \ / wn : m \ _ f \ 
\ V 0 J \ U N?M J 11 / 


(16) 


(17) 


The above motivates the following local residual operator Ri, which is defined as 
Ri : 14(ay) —* M is defined by 

Ri(v) = /(M) (l) 0^,M + {Ah) ( A'E.MUN,M — fh'j 

and the local residual norm, ||-R;|| is defined by 


non l^(«)l 

Lrt d = sup — . , 

vev h (u h) y/v T (A h AlY iS >v 

where (AhA^)^ and A^ are local sub-matrices of Ah AJ) and Ah which correspond to the 
coarse grid subdomain co l . Next, we use the local residual to construct the local test basis, 
4>on G Vh(uJi) such that 


v T (A h Al)(' ) 4 = Ri(v), Vv e V h (a;*). 

In [13], we show that if online basis functions are constructed using the second eigenvalue 
problem, then the error will decrease at a rate (1 — min K A4 n ), where A® in is the min¬ 
imum of the eigenvalues of the spectral problem defined on Wg nap (E) corresponding to 
eigenfunctions not chosen as basis. Consequently, using online basis functions, we can 
achieve the discrete inf-sup stability in one iteration provided inin B AC n > 0. 


3.5.1. Online test basis enrichment algorithm 

First, we choose an offline trial space, V 0 g and an initial offline test space, , by 
fixing the number of basis functions for each coarse neighborhood. Next, we construct a 
sequence of online test spaces and compute the multiscale solution (wm" by 

solving equation (17). The test space is constructed iteratively for m — 1, 2, 3,..., by the 
following algorithm: 
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Step 1: Find the multiscale solution in the current space. Solve for (w^ ■>,«£?) € w£>xv oS 
such that 


Step 2: For each coarse region a;*, compute the online basis, 4>ol G I4(^j) such that 

v T (A h A T h fU% = Ri(v), Vi> e V,M). 



Step 3: Enrich the test space by setting 



We remark that in each iteration, we perform the above procedure on non-overlapping 
coarse neighborhoods, see na. 

4. Numerical Results 

In this section, we present representative numerical examples. In all our examples, {xi} 
is a multiscale partition of unity. In each coarse space, we compare the Z 2 projection error 
and the Z 2 error for the multiscale solution. For simplicity, we refer to “the multiscale error” 
as the error between the multiscale solution and the exact solution, and “the projection 
error” as the error between the exact solution and its / 2 projection onto the span of the 
coarse trial space. We also assume n is a constant and b is a multiscale held. In particular, 
the velocity fields contain oscillations and cells (eddies, separatrices and/or layers) within 
a single coarse block of the discretization and, thus, we do not have a single streamline 
direction per coarse block. Fully resolved velocity solutions are shown in Figures [4] to [6] 
and in Figure [HJ The method can easily handle multiscale diffusion coefficients. The 
fine-grid problem is always chosen such that the local Peclet number is about 1 ensuring 
a stable fine discretization. All coarse discretizations have a Peclet number at least an 
order of magnitude larger than 1. We analyze the performance of the trial and test spaces 
proposed in the previous section. We pay special attention to the effect of eigenvalue 
problem on the performance of the discrete system and discuss this for each example. 


Example 1 


These first numerical examples are defined by the following diffusion and convection 
coefficients and right-hand side, 


1 


” 100 ’ 


b = a 



/ — 1 - 
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0 0.2 0.4 0.6 0.8 1 0 0.2 0.4 0.6 0.8 1 


(a) a = 2 (b) a = 4 

Figure 4: Reference solutions for Example 1. 


T^basis L 2 (projection error) 

(trial, test) Eigenproblem 1 Eigenproblem 2 



a = 2 

a = 4 

a = 2 

a = 4 

(1.1) 

8.56% 

7.06% 

11.94% 

9.60% 

(1.3) 

3.22% 

4.96% 

4.74% 

4.48% 

(1.5) 

2.85% 

4.74% 

2.90% 

5.02% 

(1.7) 

2.85%(2.85%) 

3.64%(3.52%) 

2.85%(2.85%) 

3.55%(3.52%) 

(3.1) 

9.00% 

7.58% 

11.95% 

8 .86% 

(3.3) 

3.12% 

5.22% 

5.01% 

3.96% 

(3,5) 

2.61% 

3.96% 

2.70% 

4.83% 

(3.7) 

2.60%(2.60%) 

3.41%(3.21%) 

2.61%(2.60%) 

3.25%(3.21%) 

(5.1) 

8.65% 

7.88% 

12.80% 

9.08% 

(5.3) 

2.72% 

4.97% 

4.69% 

3.35% 

(5.5) 

2.31% 

3.62% 

2.37% 

3.99% 

(5.7) 

2.31%(2.31%) 

2.89%(2.77%) 

2.31%(2.31%) 

2.79% (2.77%) 


Table 1: Errors for test space derived using Eigenproblems 1 and 2 for Example 1. Coarse and fine mesh 
sizes are H = 1/10 and h = 1/200, respectively. The projection errors are in parentheses. 


This velocity field has a cellular structure with several eddies and separatices. In the 
simulations a takes values of 2 and 4. Figure [4] depicts well-resolved fine-scale solutions 
for the chosen values of a. In both cases, we take the coarse mesh size to be H = 1/10, 
while the fine mesh size to be h = 1/200. 

Table [T| shows the impact of increasing the number of coarse basis functions as well 
as how the system converges as we increase the number of test functions included per 
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#basis test 

min{A 

Li+l} 


a = 2 

a = 4 

1 

0.3445 

0.3693 

3 

0.7273 

0.7707 

5 

0.9542 

0.9514 

7 

0.9908 

0.9919 


Table 2: Minimum eigenvalue for the test space constructed using the Eigenproblem 2 (minimal energy 
test functions) for Example 1. Coarse and fine mesh sizes are H = 1/10 and h = 1/200, respectively. 


coarse block edge. The table shows the evolution of the multiscale error as we increase 
the number of test functions for different numbers of trial functions in each coarse block. 
Each column is labeled by its corresponding value of a. Table [T| shows the performance 
of the reduced-dimensional test space constructed using the first and second eigenvalue 
problems we describe in Section 3.2.2 This table shows that 7 test functions per edge of 
the coarse mesh are enough to deliver similar multiscale and projection errors irrespective 
of a (i.e., coarse scale Peclet number) and the number of coarse basis functions used in 
each coarse block. In fact, these errors are similar even when the number of test functions 
is 5. Table [2] shows the evolution of the minimum eigenvalue for the test space constructed 
using the Eigenproblem 2 (minimal energy test functions) of Section 3.2.2 As it follows 


from the theory, for a rich enough test space with a sufficient number of multiscale test 
functions, multiscale and projection errors converge. The eigenvalue behavior shown in 
Tableland the convergence shown in Table [l] verifies that when the minimum eigenvalue 
is close to 1, the multiscale error converges to the projection error. 




(a) a = 2 


(b) a = 4 


Figure 5: Reference solutions for Example 2. 
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Example 2 

We consider the following right-hand side, and diffusion and convection coefficients: 


k = 1 / 100 , 

^ _ f + sin(187nr) cos(187 ry) + 5 cos(18\/27rx) sin(18\/27n/) 
\ — cos(187nr) sin(187n/ — 5 sin(18v^7r^) sin(18\/27n/) 

/ — 1 - 


J 


This velocity held has a cellular structure with eddies and channels, as can be seen in the 
fine-scale solutions depicted in Figure |5j These channels introduce global effects, a takes 
values of 2 and 4, as in our previous example. 5 takes value of \/2/4. For both cases, we 
take the coarse mesh size to be H = 1/10, while the fine mesh size to be h = 1/200. 


#basis 
(trial, test) 


(u) 

(1.3) 

(1.5) 

(1.7) 


(3,1) 

(3,3) 

(3,5) 

(3,7) 


(5,1) 

(5,3) 

(5,5) 

(5,7) 


L 2 (projection error) 

Eigenproblem 1 Eigenproblem 2 


a = 2 


a = 4 


a = 2 


a = 4 


9.49% 

3.11% 

2.90% 

2.90%(2.90%) 


8.57% 

5.04% 

4.73% 

3.87%(3.67%) 


12.27% 

4.62% 

2.95% 

2.90%(2.90%) 


9.99% 

4.22% 

4.62% 

3.70%(3.67%) 


10.04% 

2.97% 

2.65% 

2.65%(2.64%) 


9.09% 

4.50% 

3.85% 

3.56%(3.33%) 


12.35% 

5.09% 

2.75% 

2.65%(2.64%) 


9.33% 

3.73% 

4.52% 

3.38%(3.33%) 


9.22% 

2.67% 

2.36% 

2.36%(2.36%) 


9.25% 

4.18% 

3.43% 

3.18%(2.84%) 


13.37% 

4.60% 

2.42% 

2.36%(2.36%) 


9.98% 

3.26% 

3.84% 

2.88%(2.84%) 


Table 3: Errors for test space derived using Eigenproblems 1 and 2 for Example 2, with S = -\/2/4. Coarse 
and fine mesh sizes are H = 1/10 and h = 1/200, respectively. The projection errors are in parentheses. 


In Table [3j we increase the number of test functions and consider different numbers 

As before, we 


of trial functions for both eigenvalue problems described in Section 3.2.2 


observe that for 7 test functions, the multiscale and projection errors converge, while for 
5 test functions per edge, the errors are close. As we increase, the dimension of the coarse 
trial space, we observe a similar behavior. The eigenvalue behavior shows (see Table [4]) 
that when the eigenvalue is close to 1, the multiscale and projection errors converge. 
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#basis test 

min{A Li+ i} 


a = 2 

a = 4 

1 

0.3551 

0.2967 

3 

0.7289 

0.6761 

5 

0.9510 

0.9331 

7 

0.9909 

0.9813 


Table 4: Minimum eigenvalue for the test space constructed using the Eigenproblem 2 (minimal energy 
test functions) for Example 2 for 5 = %/2/4. Coarse and fine mesh sizes are H = 1/10 and h = 1/200, 
respectively. 



(a) a = 1/1000 (b) a = 1/2000 

Figure 6: Reference solutions for Example 3. 


Example 3 

We consider the following diffusion and convection coefficients and right-hand side. 

k — a, 



f = 1, 


where 

H = (sin(57nc) sin(67n/)/(607r)) + 0.005(a; + y). 

This velocity held again has a cellular structure with eddies and channels, as the fine- 
scale solutions show in Figure [6j In this example, ct is a diffusion coefficient and we 
take a = 1/1000 and a = 1/2000. For both cases, we take the coarse mesh size to be 
H — 1/10, while the fine mesh sizes are set to h = 1/400 and h = 1/800 for a = 1/1000 
and a = 1/2000, respectively. 
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T^basis 


L 2 (projection error) 


(trial, test) 

Eigenproblem 1 

Eigenproblem 2 


a = 1/1000 

a = 1/2000 

a = 1/1000 

a = 1/2000 

(U) 

9.59% 

32.36% 

19.09% 

34.55% 

(1.3) 

7.84% 

13.70% 

15.58% 

32.00% 

(1.5) 

7.83% 

12 .10% 

9.10% 

22.53% 

(1.7) 

7.83%(7.83%) 

11.85%(11.48%) 

7.95%(7.83%) 

14.70%(11.48%) 

(3.1) 

6 .01% 

16.24% 

16.70% 

37.12% 

(3.3) 

4.89% 

9.07% 

7.93% 

18.82% 

(3,5) 

4.88% 

7.98% 

5.33% 

11.90% 

(3.7) 

4.88%(4.88%) 

8.06%(7.58%) 

4.98%(4.88%) 

9.71%(7.58%) 

(5.1) 

4.70% 

14.77% 

14.00% 

31.46% 

(5.3) 

3.74% 

6.61% 

4.79% 

13.66% 

(5.5) 

3.72% 

6.23% 

3.80% 

7.55% 

(5.7) 

3.72% (3.72%) 

6.21%(6.15%) 

3.74%(3.72%) 

6.53%(6.15%) 


Table 5: Errors for test space derived using Eigenproblems 1 and 2 for Example 3. Discrete parameters 
used are: a = 1/1000, H = 1/10, h = 1/400 and a = 1/2000, H = 1/10, h = 1/800. The projection 
errors are in parentheses. 


#basis test 

min{A 

Li+1 } 


a = 1/1000 

a = 1/2000 

1 

0.3547 

0.3068 

3 

0.7497 

0.6304 

5 

0.9546 

0.8718 

7 

0.9952 

0.9754 


Table 6: Minimum eigenvalue for test space derived using Eigenproblem 2 for Example 3. Discrete 
parameters used are: a = 1/1000, H = 1/10, h = 1/400 and a = 1/2000, H = 1/10, h = 1/800. 

Tables [5] and [6] show a similar behavior to the one discussed in the previous two 
examples. That is, table [5] shows that if only one test function is chosen for a = 1/2000, 
the error is about 35% (when the number of test functions coarse edge is 1 and the number 
of trial functions per coarse block is 1). These errors rapidly drop to about the projection 
error as we increase the dimension of the test space. Similar behavior is observed for both 
eigen-constructions of the test space as we refine the coarse trial space. 
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#basis (trial, test) 

L 2 (projection error) 

(i.i) 

20 .12% 

(1.3) 

11.25% 

(1.5) 

4.03% 

(1.7) 

3.93%(3.93%) 

(3,1) 

19.99% 

(3.3) 

13.02% 

(3.5) 

3.31% 

(3.7) 

3.23%(3.23%) 

(5,1) 

13.93% 

(5.3) 

9.14% 

(5.5) 

2.90% 

(5.7) 

2.74%(2.70%) 


Table 7: Errors for test space derived using Eigenproblem 2 for Example 4. Coarse and fine mesh sizes 
are H = 1/10 and h = 1/200, respectively. The projection errors are in parentheses. 


Example 4 

As we remove the eddies from the flow held and make the how more channelized, the 
multiscale error grows. To expose this behavior, we take the velocity held to be 


K = 1 , 

b = 200 ^ sin ( 18 ^ 7ry ) 

/ = 1 , 

which corresponds to solving the how equations with a channelized permeability held. 
The numerical results are presented in Table [T] (the mesh sizes for the coarse and hne 
spaces are H = 1/10 and h = 1/200). We observe that the multiscale error is 20.12% for 
one trial function per coarse block and one test function per coarse interface, while the 
error reduces to the projection error of 3.93% when we select 7 test functions interface. As 
before, for 5 trial functions per coarse block, it takes 7 test functions per coarse interface to 
reduce the error to the projection error from 13.93%. For this discrete problem setup, the 
smallest eigenvalue is 0.9952 for 7 test functions per edge when minimal energy functions 


are used (Eigenproblem 2 in Section 3.2.2). 


Example 5 

We consider the following diffusion and convection coefficients, and right-hand side. 


k — a, 
b = nSJp 

f = 1, 
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Figure 7: Permeability field used to compute the transport velocity field in Example 5. The black region 
corresponds to the permeability 500 and the white region corresponds to the permeability 1. 



(a) a = 1/250 (b) a = 1/500 

Figure 8: Reference solutions for Example 5. 

where the velocity field solves this flow equation 

—V • (nVp) = 0 
Plan = xy 

Figure [7] shows the permeability held used in the above equation. The resulting velocity 
held contains channels with variable velocity in each coarse region. Figure [8] shows the 
hne-scale structure of the fully resolved velocity held. In this case a is a diffusion coefficient 
and takes values 1/250 and 1/500. In both cases, we take the coarse mesh size is set to 
H = 1/10, while the hne mesh sizes are h = 1/200 and h = 1/400 for a = 1/250 and 
a = 1/500, respectively. 
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T^basis 
(trial, test) 

L 2 (projection error) 

1 Eigenproblem 1 Eigenproblem 2 


a = 1/250 

a = 1/500 

a = 1/250 

a = 1/500 

(U) 

2 .11% 

4.64% 

2.24% 

4.56% 

(1.3) 

2.08% 

4.59% 

2.09% 

4.24% 

(1.5) 

2.07% 

4.45% 

2.07% 

4.15% 

(1.7) 

2.07%(2.07%) 

4.23%(4.04%) 

2.07%(2.07%) 

4.19%(4.04%) 

(3.1) 

1 .01% 

1.98% 

1 .68% 

3.57% 

(3.3) 

0.99% 

1.93% 

1.07% 

2.36% 

(3,5) 

0.99% 

1.93% 

1 .00% 

2.05% 

(3.7) 

0.99%(0.99%) 

1.93%(1.91%) 

0.99%(0.99%) 

2.01%(1.91%) 

(5.1) 

0.85% 

1.70% 

1.64% 

4.12% 

(5.3) 

0.75% 

1.44% 

0.84% 

1.91% 

(5.5) 

0.75% 

1.44% 

0.76% 

1.53% 

(5.7) 

0.75%(0.75%) 

1.44%(1.42%) 

0.75%(0.75%) 

1.49%(1.42%) 


Table 8: Errors for test space derived using Eigenproblems 1 amd 2 for Example 5. Discrete parameters 
used are: a = 1/250, 7? = 1/10, h = 1/200 and a = 1/500, 7? = 1/10, h = 1/400. The projection errors 
are in parentheses. 

Tables [8] and [9] show a similar behavior to that observed in the prior examples. That 
is, the multiscale error converges to the projection error as we increase the number of 
test functions per coarse edge for either eigenvalue problem and for any number of coarse 
functions in each coarse block. 


#basis test 

min{A 

Li+1 } 


a = 1/1000 

a = 1/2000 

1 

0.4106 

0.3544 

3 

0.8592 

0.7583 

5 

0.9828 

0.9535 

7 

0.9985 

0.9919 


Table 9: Minimum eigenvalue for test space derived using Eigenproblem 2 for Example 3. Discrete 
parameters used are:a = 1/250,7? = 1/10, h = 1/200 and a = 1/500,7? = 1/10, h = 1/400 


4-1. Numerical Result for online test basis enrichment 

In this section, we present some numerical results, which use online test basis functions 
to stabilize the system. In Table 10, we show the convergence history for the online test 
basis enrichment for Example 1, while in Table 11 we show the convergence history of the 
online test basis enrichment for Example 4. In these two cases, with only one iteration, 
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#basis 
(trial, test) 

#iter 

Eigenproblem 1 
a = 2 a = 4 

Eigenproblem 2 
a = 2 a = 4 


0 

8.56% 

7.06% 

11.94% 

9.60% 

(u) 

1 

2.89% 

3.79% 

2.96% 

4.30% 


2 

2.85%(2.85%) 

3.52%(3.52%) 

2.85%(2.85%) 

3.52%(3.52%) 


0 

3.22% 

4.96% 

4.74% 

4.48% 

(1.3) 

1 

2.85% 

3.54% 

2 .86% 

3.63% 


2 

2.85%(2.85%) 

3.52%(3.52%) 

2.85%(2.85%) 

3.52%(3.52%) 


0 

8.65% 

7.88% 

12.80% 

9.08% 

(5.1) 

1 

2.33% 

2.97% 

2.58% 

3.29% 


2 

2.31%(2.31%) 

2.77% (2.77%) 

2.32%(2.31%) 

2.78%(2.77%) 


0 

2.72% 

4.97% 

4.69% 

3.35% 

(5,3) 

1 

2.31% 

2.79% 

2.33% 

2.81% 


2 

2.31%(2.31%) 

2.77% (2.77%) 

2.31%(2.31%) 

2.77%(2.77%) 


Table 10: Error evolution as online basis functions are added to the system (test space derived using 
Eigenproblems 1 and 2 for Example 1). Coarse and fine mesh sizes are H = 1/10 and h = 1/200, 
respectively. The projection errors are in parentheses. 

the multiscale error becomes similar to the projection error. In the second iteration, the 
multiscale error converges to the projection error. 


#basis 

#iter 


(trial, test) 




0 

20 .12% 

(u) 

1 

3.93% 


2 

3.92%(3.92%) 


0 

11.15% 

(1,3) 

1 

3.92% 


2 

3.92%(3.92%) 


0 

13.93% 

(5,1) 

1 

3.24% 


2 

2.72% (2.70%) 


0 

9.14% 

(5,3) 

1 

2.74% 


2 

2.70% (2.70%) 


Table 11: Error evolution as online basis functions are added to the system (test space derived using 
Eigenproblem 2 for Example 4). Coarse and fine mesh sizes are H = 1/10 and h = 1/200, respectively. 
The projection errors are in parentheses. 
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5. Conclusions 

In this paper, we study multiscale methods for convection-dominated diffusion with 
heterogeneous convective velocity fields. This stabilization generalizes the approaches 
described in B3I to multiscale problems. To construct this stabilization we reformulate 
overall problem in mixed form. The auxiliary variable we introduce plays the role of the 
test function. We describe the multiscale spaces we use for the test and trial spaces, 
which are built using GMsFEM framework. First, we construct snapshots spaces. For 
the test variable, we propose local snapshot spaces. Furthermore, we propose a local 
spectral decomposition following our recent work da, where we consider minimum energy 
snapshot functions. We discuss the discrete stability of the system and its relation to 
the approximation properties of the velocity field. The resulting approximation error is 
minimized within our multiscale framework by selecting a few multiscale basis functions. 
Our numerical results show that we can stabilize the system using a few test functions 
for a given trial space. We describe and analyze several relevant numerical examples that 
validate our theoretical results. 
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